#Tai-Hsien Ou Yang 06092012
FILE_LOCATION   ="/home/to2232/dataset/"
setwd(FILE_LOCATION)
pr<-c()
bk<-c()

parchild<-c()
bkg<-c()
for(k in 1:22)
{
	#"K0003" "K5826" "792429" "245377441" "271.3" "0" "0" "1" "1"
	FILE_OUT_PR           = paste("p.",k,".txt", sep="")
	FILE_OUT_NR           = paste("n.",k,".txt", sep="")
	if (file.exists(FILE_OUT_PR)){
	#load(FILE_OUT_PR)
	parchild<-scan(FILE_OUT_PR , list(ch="",pr="",sgs=0,sge=0,len=0,hom1=0,hom2=0,note1,note2))
	pr<-rbind(pr, parchild)
	}
	if (file.exists(FILE_OUT_NR)){
	#load(FILE_OUT_NR)
	bkg<-scan(FILE_OUT_PR , list(ch="",pr="",sgs=0,sge=0,len=0,hom1=0,hom2=0,note1,note2))
	bk<-rbind(bk, bkg)
	}
	
	
}


cat("Draw Histogram\n")

scale.x=1000
scale.y=2000



#Total Length Histogram
pos.len<-hist(pr$len,seq(0, scale.x, 1))
jpeg('hist_lenp.jpg')
plot(pos.len)
dev.off()

neg.len<-hist(bk$len,seq(0, scale.x, 1))
jpeg('hist_lenn.jpg')
plot(neg.len)
dev.off()

tiff('hist_len.tiff',width = scale.y, height = scale.y, units = "px")
plot( pos.len, col=rgb(0,0,1,1/4), xlim=c(0,scale.x),main="IBD Length",sub = "OuYang Apr 20",xlab="Length (cM)")  # first histogram
plot( neg.len, col=rgb(1,0,0,1/4), xlim=c(0,scale.x),, add=T)  # second
#plot( rslt_2, col=rgb(0,1,0,1/4), xlim=c(0,scale.x),ylim=c(0,1.5e6), add=T)  # second
dev.off()

pos.nm<-hist(pr$sge-pr$sgs,seq(0, scale.x, 10))
jpeg('hist_mnp.jpg')
plot(pos.nm)
dev.off()

#Segment Histogram
neg.mn<-hist(bkg$sge-bkg$sgs,seq(0, scale.x, 10))
jpeg('hist_mnn.jpg')
plot(neg.mn)
dev.off()
tiff('hist_mn.tiff',width = scale.y, height = scale.y, units = "px")
plot( pos.nm, col=rgb(0,0,1,1/4), xlim=c(0,scale.y),main="IBD Segment",sub = "OuYang Apr 20",xlab="Position")  # first histogram
plot( neg.mn, col=rgb(1,0,0,1/4), xlim=c(0,scale.y),, add=T)  # second
#plot( rslt_2, col=rgb(0,1,0,1/4), xlim=c(0,scale.y),ylim=c(0,1.5e6), add=T)  # second
dev.off()
